home *** CD-ROM | disk | FTP | other *** search
/ APDL Eductation Resources / APDL Eductation Resources.iso / programs / electronic / rlab / !RLaB / examples / curve_fit < prev    next >
Encoding:
Text File  |  1995-01-24  |  1.1 KB  |  64 lines
  1. //
  2. // curve_fit.r
  3. //
  4.  
  5. // Use a minimization scheme to fit a model to a set of data
  6. // points.
  7.  
  8. // Create the data...
  9.  
  10. require rem
  11.  
  12. x = (0:10:.1)';
  13. y = 3.2.*1.5.^x;
  14.  
  15. ptitle ("Original Contrived Data");
  16. plot ([x,y]);
  17. pause ();
  18.  
  19. // Now add some random-ness
  20.  
  21. rand ("uniform", -max(abs(y))/5, max(abs(y))/5);
  22. y_orig = y;
  23. y = y + rand (size (y));
  24.  
  25. ptitle ("Original Data with some Random-ness Added");
  26. plot ([x,y]);
  27. pause ();
  28.  
  29. // Now make up a model function
  30.  
  31. fmodel = function (a, b, x) { return a.*b.^x; }
  32.  
  33. // Now make up a cost function
  34.  
  35. fcost = function ( A ) 
  36. {
  37.   global (x, y)
  38.   return -sum ( abs ((y - fmodel (A[1], A[2], x))).*x); 
  39. }
  40.  
  41. // Now load the min/max - imizer
  42.  
  43. rfile mdsmax
  44. rfile nmsmax
  45.  
  46. #fit = nmsmax (fcost, [1, 1]);
  47. fit = mdsmax (fcost, [1, 1]);
  48.  
  49. // Plot the result
  50.  
  51. ptitle ("Original+Random-ness and Fitted Model");
  52. plot ([ x, y, fmodel (fit.x[1], fit.x[2], x) ]);
  53. pause ();
  54.  
  55. ptitle ("Original and Fitted Model");
  56. plot ([ x, y_orig, fmodel (fit.x[1], fit.x[2], x) ]);
  57.  
  58. "\tNow compare original model coefficients and fitted"
  59.  
  60. original = [3.2, 1.5]
  61. fitted = fit.x
  62.  
  63. difference = [3.2, 1.5] - fit.x
  64.